Wykonali: Majka Miezianko, Marek Czerwonka, Antoni Wadowski

Wykrywanie obecności choroby serca na podstawie wybranych czynników

Podstawy teoretyczne tworzenia sieci neuronowej

Wstęp

Przeznaczeniem naszego modelu jest odwzorowywanie wartości opisywanej (w tym przypadku występowanie choroby serca) na podstawie danych wejściowych pochodzących ze zbioru danych o wspólnych aspektach (np. parametry związane z kondycją zdrowotną człowieka) . Jeśli zdefiniujemy model jako funkcję przyjmującą na wejściu parametry, możemy tak go dopasować, by w miarę możliwości jak najbardziej optymalnie opisywał dane. W uogólnieniu służy do tego następująca procedura:

  1. Wielokrotne przekazywanie obserwacji do modelu i zapisywanie wartości obliczanych w kroku w przód (forward propagation)
  2. Obliczanie wartości straty (loss), reprezentującej, jak odległe są predykcje modelu od oczekiwanych danych wyjściowych
  3. Używanie wartości obliczonych w forward propagation i opartych na regule łańcuchowej obliczeń, aby ustalić w jakim stopniu każdy z wejściowych parametrów wpływa na wartość straty.
  4. Modyfikowanie wartości parametrów, tak aby wartość straty zmalała po przekazaniu do modelu następnego zbioru obserwacji.

Najważniejszy dla nas jest czwarty krok. Chcemy mieć dowolną liczbę warstw, dowolnego typu. Jeśli jednak zmodyfikujemy / dodamy / usuniemy jedną warstwę z sieci, wyjście sieci również ulegnie zmianie, co spowoduje zmianę błędu, który zmieni pochodną błędu w odniesieniu do parametrów. Musimy być w stanie obliczyć pochodne niezależnie od architektury sieci, niezależnie od funkcji aktywacji, niezależnie od strat, z których korzystamy. Aby to osiągnąć, zaimplementowaliśmy w projekcie każdą warstwę osobno.

Fully connected layers: warstwy w pełni połączone

W sieciach tego typu każda jednostka wyjściowa jest wyliczana jako suma ważona wszystkich wejść. Określenie "w pełni połączona" oddaje sposób działania sieci: każde wyjście jest połączone ze wszystkimi wejściami. Większość publikacji przedstawia sieci w pełni połączone za pomocą zapisu macierzowego. W tym przypadku mnożymy wektor wejść przez macierz wag W i uzyskujemy w ten sposób wektory wyjść:

y = Wx

Ponieważ mnożenie macierzy jest operacją liniową, sieć zawierająca jedynie mnożenie macierzy byłaby w stanie uczyć się tylko odwzorowań liniowych. Aby sieci neuronowe miały większe możliwości wyrazu, mnożenie macierzy jest uzupełniane przez nieliniową funkcję aktywacji. Może to być dowolna funkcja różniczkowalna, jednak powszechnie stosuje się tylko kilka takich funkcji. Do niedawna najczęściej stosowaną funkcją aktywacji był tangens hiperboliczny, który również został zaimplementowany w jednej z klas naszego projektu. Problem związany ze stosowaniem tej funkcji polega na tym, że staje się ona bardzo "płaska" dla wartości wejściowych znacząco odbiegających od zera. W efekcie gradient tej funkcji jest bardzo mały, co oznacza, że uzyskanie zmiany zachowania sieci może zająć bardzo dużo czasu.

Warstwy umieszczone wewnątrz sieci są nazywane warstwami ukrytymi, a poszczególne wyjścia używane w warstwach ukrytych - jednostkami ukrytymi. Określenie "ukryty" pochodzi stąd, że jednostki te nie są widoczne bezpośrednio na zewnątrz sieci, tak jak jej wejścia i wyjścia.

Forward propagation: krok w przód

Na początku należy podkreślić jedną ważną kwestię, na której oparte są modele deep learningowe: dane wyjściowe jednej warstwy są danymi wejściowymi kolejnej warstwy. Nazywa się to propagacją do przodu. Zasadniczo przekazujemy dane wejściowe do pierwszej warstwy, a następnie dane wyjściowe każdej warstwy stają się danymi wejściowymi kolejnej warstwy, aż do końca sieci.

W praktyce wygląda to tak, że wejście neuronu jest równe sumie ważonych wyników ze wszystkich neuronów w poprzedniej warstwie. Każde wejście jest mnożone przez wagę związaną z synapsą łączącą wejście z bieżącym neuronem. Jeśli w poprzedniej warstwie są 3 neurony, każdy neuron w bieżącej warstwie będzie miał 3 różne wagi.

Następnie porównując wynik sieci (Y empiryczne) z żądanym wyjściem (powiedzmy Y teoretyczne) możemy obliczyć błąd E. Celem jest zminimalizowanie tego błędu poprzez zmianę parametrów w sieci. Ten proces nazywany jest back propagation (propagacja wsteczna).

Wartość każdego neuronu wyjściowego można obliczyć w następujący sposób:

Mając macierze, możemy obliczyć tę formułę dla każdego neuronu wyjściowego:

Gradient descent

Gradient descent jest algorytmem optymalizacji używanym do minimalizowania niektórych funkcji poprzez iteracyjne poruszanie się w kierunku najbardziej stromego spadku (ujemny gradient). W uczeniu maszynowym używamy gradient descent, aby zaktualizować parametry naszego modelu. Parametry odnoszą się do wag w sieciach neuronowych.

Podstawową intuicję związaną z gradient descent można zilustrować hipotetycznym scenariuszem. Osoba utknęła w górach i próbuje zejść (tj. próbuje znaleźć globalne minimum). Panuje gęsta mgła, przez co widoczność jest wyjątkowo niska. Dlatego ścieżka w dół nie jest widoczna, więc aby znaleźć minimum, musi skorzystać z lokalnych informacji. Potrafi skorzystać z metody gradient descent, która polega na obserwowaniu stromości wzniesienia w jej aktualnej pozycji, a następnie zejściu w kierunku najbardziej stromego zbocza (czyli w dół). Gdyby próbowała znaleźć szczyt góry (tj. maksimum), to podążałaby w kierunku najbardziej stromego wzniesienia (tj. pod górę). Korzystając z tej metody, w końcu trafiłaby w dolne partie góry lub prawdopodobnie utknęłaby w jakiejś dziurze (tj. lokalnym minimum), jak górskie jezioro. Jednak załóżmy również, że stromość wzgórza nie jest od razu oczywista przy prostej obserwacji, ale raczej wymaga wyrafinowanego instrumentu do pomiaru, który dana osoba ma w tej chwili. Pomiar nachylenia wzgórza za pomocą przyrządu zajmuje trochę czasu, dlatego powinna zminimalizować użycie przyrządu, jeśli chciałaby zejść z góry przed zachodem słońca. Trudność polega więc na wyborze częstotliwości, z jaką powinna mierzyć nachylenie wzgórza, aby nie zboczyć z trasy.

W tej analogii osoba reprezentuje algorytm, a ścieżka w dół góry reprezentuje sekwencję ustawień parametrów, które będzie badał algorytm. Strome wzniesienie reprezentuje nachylenie funkcji błędu w tym punkcie. Instrumentem używanym do pomiaru stromości jest różniczkowanie (nachylenie powierzchni błędu można obliczyć, biorąc pochodną funkcji błędu). Kierunek, w którym zdecyduje się podróżować, jest zgodny z gradientem błędu w tym punkcie. Czas, jaki mija przed wykonaniem kolejnego pomiaru, to rozmiar kroku - learning rate. O learning rate i metodzie jego optymalizacji wspomniemy jeszcze w dalszej części pracy.

Backproagation: krok wstecz

Cel propagacji wstecznej jest prosty: dostosowanie każdej wagi w sieci proporcjonalnie do tego, jak bardzo przyczynia się ona do ogólnego błędu. Jeśli iteracyjnie zredukujemy każdy błąd wagi, ostatecznie otrzymamy serię wag, które dadzą dobre prognozy.

Jak wspomnieliśmy powyżej, forward propagation można postrzegać jako długą serię zagnieżdżonych równań. Jeśli spojrzymy na forward propagation właśnie w ten sposób, to backpropagation jest po prostu zastosowaniem reguły łańcuchowej (chain rule) w celu znalezienia pochodnych błędu w odniesieniu do dowolnej zmiennej w zagnieżdżonym równaniu.

Nie zapominajmy, że dane wyjściowe jednej warstwy są danymi wejściowymi kolejnej warstwy. Co oznacza, że ∂E / ∂X dla jednej warstwy wynosi ∂E / ∂Y dla poprzedniej warstwy. Możemy znowu odwołać się do chain rule:

Krótkie omówienie matematyki oraz wzorów:

Załóżmy, że mamy macierz zawierającą pochodną błędu w odniesieniu do outputu tej warstwy (∂E / ∂Y).

Potrzebujemy zatem:

Obliczmy ∂E / ∂W. Ta macierz powinna mieć taki sam rozmiar jak W: i x j, gdzie i to liczba neuronów wejściowych, a j to liczba neuronów wyjściowych. Potrzebujemy jednego gradientu na każdą wagę:

∂E / ∂W

Korzystając z podanej wcześniej reguły łańcuchowej, możemy napisać:

Zatem:

∂E/∂B

Powyżej przedstawiliśmy pierwszy z wzorów, który został wykorzystany do aktualizacji wag. Oto kolejny z nich - dotyczący ∂E/∂B - pochodnej z błędu po biasie:

Ponownie ∂E / ∂B musi mieć taki sam rozmiar jak B - jeden gradient na każdy bias. Możemy raz jeszcze użyć reguły łańcuchowej:

Podsumowując:

∂E / ∂X

Pozostaje nam jeszcze nam ∂E / ∂X, które jest bardzo ważne, ponieważ „zachowuje się” jak ∂E / ∂Y dla warstwy wcześniejszej.

Jak we wcześniejszych przykładach i tym razem zastosujemy chain rule:

Cała macierz wygląda następująco:

Podsumowanie

W taki sposób wyprowadziliśmy niezbędne wzory dla podstawowej sieci typu fully conected.

Kolejne zastosowane przez nas formuły dla funkcji błędu, aktywacji czy optymalizatorów zostaną zaprezentowane na etapie aich analizy porównawczej.

Mini-batch gradient descent - MGD

Optymalizator to nic innego jak algorytm lub metody używane do zmiany atrybutów sieci neuronowych, takich jak wagi i szybkość uczenia się, w celu zmniejszenia strat. Innymi słowy, to, w jaki sposób powinniśmy zmieniać wagi i szybkość uczenia się naszej sieci neuronowej, aby zmniejszyć błąd, określa używany przez nas optymalizator.

Wyróżniamy różne typy optymalizatorów, a wśród nich wspomnieć można chociażby o optymalizatorach gradient descent tj.:

Wspomniane optymalizatory to po prostu inne sposoby na przekazywanie danych w sieci.

W naszym modelu zaimplementowaliśmy oraz wykorzystaliśmy pierwszy z nich i to właśnie od niego zaczniemy omówienie zastosowanych optymalizatorów. MGD to odmiana algorytmu gradient descent, który dzieli treningowy zestaw danych na małe partie, które są używane do obliczania błędu modelu i aktualizowania współczynników modelu.

Możemy to rozumieć w następujący sposób: przypuśćmy, że mamy 1000 rekordów, a rozmiar batcha = 50: wybieramy zatem losowo 50 rekordów, a następnie obliczamy sumę błędu i przekazujemy ją dalej w celu znalezienia ∂E / ∂W. Batche są tworzone w oparciu o losowy wybór zbiorów danych.

Rozmiary mini-batchy są często dostosowywane do aspektu architektury obliczeniowej, na której wykonywana jest implementacja. Na przykład potęga dwóch, która spełnia wymagania dotyczące pamięci sprzętu GPU lub procesora, takiego jak 32, 64, 128, 256 i tak dalej.

Często zalecaną wartością dla rozmiaru batchy jest 32 i to właśnie taka wartość została wykorzystana w naszym modelu. Poniżej przedstawimy kilka cytatów z artykułów naukowych, w których rekomendowano takie postępowanie:

"[batch size] is typically chosen between 1 and a few hundreds, e.g. [batch size] = 32 is a good default value, with values above 10 taking advantage of the speedup of matrix-matrix products over matrix-vector products."

"Practical recommendations for gradient-based training of deep architectures" [Yoshua Bengio, Cornell University]

"The presented results confirm that using small batch sizes achieves the best training stability and generalization performance, for a given computational cost, across a wide range of experiments. In all cases the best results have been obtained with batch sizes m = 32 or smaller.

"Revisiting Small Batch Training for Deep Neural Networks [Dominic Masters, Carlo Luschi, Cornell University]

W dalszej częsci pracy omówimy również inny wykorzystany przez nas optymalizator.

Opis podjętego problemu oraz datasetu

Współczynnik śmiertelności rośnie z dnia na dzień na całym świecie. Przyczyną tego - oprócz oczywistych pandemicznych względów - może być w dużej mierze także wzrost liczby pacjentów z chorobami układu krążenia. Zdarza się, że lekarze nie zawsze są w stanie postawić prawidłową diagnozę pacjentów zgłaszających dolegliwości ze strony klatki piersiowej. W naszej pracy zaprojektowaliśmy sieć, która potencjalnie mogłaby pomóc w skutecznej diagnostyce pacjenta oraz eliminacji błędnej diagnozy. W tym "eksperymencie" wykorzystaliśmy zestaw danych dotyczący chorób serca Heart Disease UCI. Zbiór danych obejmuje trzynaście cech, które są istotne w diagnostyce chorób serca. System jest modelowany w wielowarstwowej sieci neuronowej wyszkolonej za pomocą propagacji wstecznej i symulowany za pomocą jednokierunkowej sieci neuronowej (feedforward).

Powód wyboru tego zagadnienia

Powód wyboru tego zagadnienia częściowo został już przedstawiony przez nas powyżej. Według Światowej Organizacji Zdrowia (WHO) 17,9 mln ludzi umiera każdego roku z powodu chorób układu krążenia, co odpowiada 31% wszystkich zgonów na świecie. Wiąże się to z koniecznością posiadania niedrogiego systemu, który jest w stanie dokonać wstępnej oceny pacjenta na podstawie stosunkowo prostych testów medycznych, na które każdy może sobie pozwolić.

Opis datasetu

Tak jak wspomnieliśmy, w naszej pracy korzystamy ze zbioru danych o chorobach serca z repozytorium uczenia maszynowego. Chociaż zbiór danych ma 76 surowych atrybutów, zdecydowano, że tylko 14 z nich jest faktycznie użytecznych. Czternaste pole wskazuje na występowanie choroby serca u pacjenta.

Importowanie bibliotek

Importowanie danych

Sprawdzanie, czy występują braki danych

Statystyka opisowa

Możemy zauważyć, że ten zbiór danych zawiera 303 obserwacje i nie ma w nich braków danych. Wyróżniamy w nim łącznie 13 features (czynników objaśniających) oraz jedną kolumnę, która jest targetem, czyli zmienną zero-jedynkową. Zero oznacza brak choroby serca (138 obserwacji), zaś jeden występowanie choroby serca (165 obserwacji).

Warto również zwrócić uwagę, że skala zmiennych jest zróżnicowana — podczas gdy maksymalna wartość 'age' wynosi 77, maksymalna wartość 'chol' równa się 256.

Wizualizacja danych

W celu lepszego zrozumienia oraz objaśnienia datasetu wykonamy poniżej wizualizację danych.

Patrząc na powyższą macierz korelacji, łatwo zauważyć, że kilka cech ma ujemną korelację z wartością docelową, podczas gdy inne mają dodatnią.

Zauważamy, że najbardziej skorelowane z targetem są zmienne:

Najbardziej skorelowane pomiędzy sobą są zaś zmienne:

Jakie ma to znaczenie dla naszej sieci?

Jak można przeczytać w artykule "The Precise Effect of Multicollinearity on Classification Prediction" [Mary G. Lieberman, John D. Morris, Florida Atlantic University] współliniowość nie jest aż tak bardzo istotna dla problemów klasyfikacyjnych z użyciem walidacji krzyżowej jak na przykład w przypadku regresji. Według G. Lieberman i J. Morrisa, jeśli celem jest dokładność prognozowania, to współliniowość nie ma znaczenia. Jeśli, z drugiej strony, meritum naszych badań nie jest skupione wyłącznie na dokładności przewidywania modelu, a bardziej na znaczeniu zmiennej w modelu, to oczywiście współliniowość może utrudniać rozważania.

Warto zatem zwrócić uwagę, iż zmienne w naszym modelu nie są silnie skorelowane.

Sprawdzamy jaki jest rozkład Y - targetu wśród poszczególnych czynników

Możemy zauważyć z powyższych wykresów, że więcej kobiet niż mężczyzn zachorowało na choroby serca.

Histogramy zmiennych:

Jak możemy zauważyć, wszystkie atrybuty ciągłe (oprócz oldpeak) mają rozkład zbliżony do normalnego. W przypadku atrybutu age obserwujemy, że większość obserwacji mieści się w przedziale 47–61 lat.

W tym przypadku punkty skupiają się wokół wartości 170 i 130 odpowiednio na osi x i y.

One hot encoding

Algorytmy uczenia maszynowego nie mogą przetwarzać danych kategorialnych ani tekstowych, chyba że zostały przekonwertowane na liczby. One hot encoding odwzorowuje wartości kategorialne na wartości całkowite, które są reprezentowane jako wektor binarny, z których wszystkie są wartościami zerowymi, z wyjątkiem indeksu, który jest ustawiony na 1.

Standaryzacja danych

Normalizacja to przeskalowanie danych z pierwotnego zakresu, tak aby wszystkie wartości mieściły się w pewnym zakresie, zwykle między 0 a 1. Znormalizowane dane są niezbędne w uczeniu maszynowym.

Model właściwy sieci neuronowej

Import stworzonych przez nas klas:

Podział na zbiór treningowy i uczący

Warto na tym etapie wspomnieć również o parametrze, występującym w stworzonym przez nas modelu sieci neuronowej - proc_walidacyjny, który odnosi się oczywiście do podziału danych na podzbiór walidacyjny. Pozwolimy sobie zatem krótko nakreślić ideę jego stosowania.

Często można się również z nazwą „development set” lub „dev set”. Podkreśla ona, że zestaw danych walidacyjnych służy do testowania wyników podczas doboru modelu czy hiperparametrów na etapie rozwoju wybranego algorytmu.

Poza tym będziemy posługiwać się oczywiście, wspomnianymi powyżej: zbiorem testowym oraz treningowym.

Oto jak prezentuje się nasz ANN:

Jak możemy zauważyć accuracy stworzonego przez nas modelu wynosi 98.36%. Accuracy, czy inaczej mówiąc dokładność, to jedna z miar służących do oceny modeli klasyfikacyjnych. Nieformalnie, dokładność to procent poprawnie przewidzianych przez model prognoz. Formalnie ma ona następującą definicję:

$ Accuracy = \frac{\text{poprawne predykcje}}{\text{liczba wszystkich predykcji}} $

Błąd naszego modelu wynosi natomiast mniej niż 0.05.

Przedstawione powyżej wyniki możemy uznać za względnie satysfakcjonujące. Warto jednakże wspomnieć o tym, czemu je zawdzięczamy oraz jakie konkretnie parametry zostały przez nas wykorzystane. Poniżej prezentujemy zatem zwięzłą specyfikację sieci, a jej poszczególne parametry omówimy w części porównawczej.

Parametry ANN:

W kolejnych krokach omówimy wpływ poszczególnych parametrów oraz ich znaczenie dla wyników i dokładności naszej sieci neuronowej.

Analiza wpływu wybranych czynników na dokładność sieci neuronowej

Analiza dokładności sieci w zależności od ilości warstw

Chociaż liczba ukrytych warstw i liczba neuronów na warstwę są kluczowe dla modeli uczenia maszynowego, nie ma ogólnej zasady ich określania w celu uzyskania najlepszej wydajności. Zwykle jest to metoda heurystyczna i może mieć na nią wpływ wiele czynników, takich jak liczba danych treningowych i zastosowane funkcje.

Co prawda trudno powiedzieć, aby naszą metodę doboru można było nazwać heurystyczną - jednak z pewnością bazuje ona na empirycznych spostrzeżeniach.

W naszym przypadku najlepiej sprawdziła się liczba warstw równa 3.

Analiza dokładności sieci w zależności od występowania optymalizatora learning rate

Jak wspominaliśmy na początku, mini-batch gradient descent to nie jedyny optymalizator, który został zastosowany w naszym modelu. Kolejnym z nich jest optymalizator learning rate, znany jako Adaptive Moment Estimation (Adam), którego ideę postaramy się przedstawić poniżej.

Adam został po raz pierwszy zaprezentowany przez Diederika Kingmę z OpenAI i Jimmy'ego Ba z University of Toronto w artykule ICLR z 2015 roku zatytułowanym „Adam: A Method for Stochastic Optimization”.

Przedstawiając algorytm, autorzy wymienili wiele korzyści, wynikających z jego zastosowania - my przytoczymy jedynie parę z nich:

Adam różni się od klasycznego SGD. Stochastic Gradient Descent utrzymuje jeden współczynnik uczenia się (learning rate) dla wszystkich aktualizacji wagi i nie zmienia się on podczas treningu. W Adamie natomiast obliczane są indywidualne współczynniki adaptacyjnego uczenia się dla różnych parametrów. Jego nazwa wywodzi się od Adaptive Moment Estimation, a powód, dla którego tak się nazywa, jest taki, że Adam używa oszacowań pierwszego i drugiego momentu gradientu, aby dostosować szybkość uczenia się dla każdej wagi sieci neuronowej. Co to jest moment? N-ty moment zmiennej losowej definiuje się jako wartość oczekiwaną tej zmiennej do potęgi n. Bardziej formalnie:

Należy zauważyć, że gradient funkcji kosztu sieci neuronowej można uznać za zmienną losową, ponieważ zwykle ocenia się go na małej, losowej porcji danych. Pierwszy moment to średnia, a drugi to tzw. uncentered variance - co oznacza, że nie odejmujemy średniej podczas obliczania wariancji. W dalszej części pokażemy, jak wykorzystujemy te wartości, teraz musimy jednak zdecydować, jak je w ogóle zdobyć.

Aby oszacować momenty, Adam wykorzystuje wykładnicze średnie ruchome, obliczone na gradiencie oszacowanym na bieżącym mini-batchu:

Gdzie m i v są średnimi ruchomymi, g jest gradientem na bieżącym mini-batchu, a beta - nowo wprowadzonymi hiperparametrami algorytmu. Parametry beta mają naprawdę dobre wartości domyślne wynoszące odpowiednio 0,9 i 0,999 i z tego, co udało nam się dowiedzieć - prawie nikt nigdy nie zmienia tych wartości. Wektory średnich ruchomych są inicjowane zerami w pierwszej iteracji.

Aby zobaczyć, jak te wartości korelują z momentem zdefiniowanym w pierwszym równaniu, przyjrzyjmy się wartościom oczekiwanym naszych średnich ruchomych. Ponieważ m i v są estymatorami pierwszego i drugiego momentu, chcemy uzyskać następującą własność:

Wartości oczekiwane estymatorów powinny być równe parametrowi, który próbujemy oszacować - a tak się składa, że parametr w naszym przypadku jest również wartością oczekiwaną. Gdyby te własności były prawdziwe, oznaczałoby to, że mamy estymatory nieobciążone. Jak jednak zaraz pokażemy, nie jest to prawda dla średnich ruchomych. Ponieważ wartości średnie inicjalizujemy zerami, estymator jest obciążony i jego obciążenie dąży do zera. Udowodnimy to dla m (dowód na v byłby analogiczny), obliczając m dla pierwszych gradientów. Poniżej rozwiniemy kilka wartości m, aby zobaczyć regułę, którą wykorzystamy w dalszym kroku:

Jak widać, im „dalej” zwiększamy wartość m, w tym mniejszym stopniu pierwsze wartości gradientów wpływają na ogólną wartość, ponieważ są one mnożone przez coraz mniejsze współczynniki beta. Zauważając tę zależność, możemy zapisać wzór na naszą średnią ruchomą:

Teraz przyjrzymy się wartości oczekiwanej m, abyśmy mogli zobaczyć, jak ma się ona do prawdziwego pierwszego momentu i skorygować ich rozbieżność:

Jak możemy zauważyć, w wyniku przekształceń matematycznych we wzorze pojawił się błąd C. W tym równaniu powinniśmy zatem zwrócić uwagę na dwie rzeczy:

Jedyne, co pozostaje do zrobienia, to wykorzystanie tych średnich ruchomych do skalowania learning rate indywidualnie dla każdego parametru. Sposób, w jaki to się robi w Adamie jest bardzo prosty. Aby zaktualizować wagę, wykonujemy następujące czynności:

Gdzie w to wagi modelu, eta (wygląda jak litera n) to rozmiar kroku. I to wszystko. Teraz skupimy się na meritum, do którego zmierzać miały nasze wyjaśnienia: czy nasz model jest "beneficjentem" Adama?

Na poniższych wykresach zestawimy ze sobą dwa takie same modele, które różni wyłącznie obecność optymalizatora.

Jak możemy zauważyć sieć, w której został zastosowany Adam wypada znacznie lepiej niż ta bez optymalizatora. Różnice są znaczne, gdyż w przypadku dokładności treningowej mamy do czynienia z około 98% w przypadku Adama i niecałym 70% bez zastosowania optymalizatora.

Analiza dokładności sieci w zależności od ilość neuronów

Liczba neuronów użytych w warstwach ukrytych została uzyskana poprzez eksperymenty, mające na celu określenie jak najlepszej zdolności predykcyjnej modelu.

Jak możemy zauważyć mała liczba neuronów nie najlepiej odwzorowuje zależności w naszym modelu. Może to wynikać z faktu, że wykorzystany przez nas surowy dataset posiada 13 cech (25 zmiennych po one hot encoding) i każda z nich powinna mieć (naszym zdaniem) pełną reprezentację neuronów. Empirycznie najlepiej sprawdziła się ilość neuronów równa 64 lub 128. W naszym modelu postanowiliśmy wykorzystać zatem 64 neurony w warstwie ukrytej.

Analiza dokładności sieci w zależności od funkcji błędu

Zanim zaczniemy porównanie wyników z wykorzystaniem poszczególnych funkcji błędu - Cross Entropy oraz MSE - postaramy się poniżej przedstawić nieco więcej informacji na temat pierwszej z nich.

Binary crossentropy to funkcja błędu, która jest (z dużym powodzeniem) używana w zadaniach klasyfikacji binarnej. Są to zadania odpowiadające na pytanie z tylko dwoma możliwościami (tak lub nie, A lub B, 0 lub 1, lewo lub prawo czy jak w naszym przypadku: chory lub zdrowy).

Cross entropy oblicza błąd danego przykładu poprzez następującą średnią:

gdzie y z daszkiem to i-ty skalar outputu, y to odpowiadająca mu wartość targetu, a output size to liczba wartości skalarnych w outpucie modelu.

Warto dodać, że zazwyczaj stosowanymi funkcjami aktywacji przy zastosowaniu entropii krzyżowej do obliczania błędu są: sigmoid (w naszym przypadku zdecydowaliśmy się na użycie właśnie tej funkcji) lub softmax - ponieważ cross entropy mierzy wydajność modelu klasyfikacyjnego, którego wynikiem jest wartość prawdopodobieństwa z zakresu od 0 do 1. Błąd entropii krzyżowej rośnie, gdy przewidywane prawdopodobieństwo odbiega od rzeczywistej etykiety. Zatem przewidywanie prawdopodobieństwa 0,012, gdy rzeczywista etykieta obserwacji wynosi 1 skutkowałoby wysoką wartością straty.

Powyższy wykres przedstawia zakres możliwych wartości błędu przy prawdziwej obserwacji (np. CzyZdrowy = 1). Gdy przewidywane prawdopodobieństwo zbliża się do 1, błąd powoli maleje. Jednak wraz ze spadkiem przewidywanego prawdopodobieństwa strata logarytmu gwałtownie rośnie. Log loss "kara" zatem szczególnie te predykcje, które są błędne.

Co więcej, z publikacji naukowych wynika, że sigmoid + mse zbiega wolniej w porównaniu z sigmoidem + entropią krzyżową ze względu na problem z zanikaniem gradientu. Przeczytać o tym można między innymi w publikacji "Cross-Entropy vs. Squared Error Training: a Theoretical and Experimental Comparison" [Pavel Golik, Patrick Doetsch, Hermann Ney,Aachen University].

Poniżej sprawdzimy jak te teoretyczne założenia sprawdzają się na naszym zagadnieniu klasyfikacyjnym.

Jak możemy zauważyć nasze przypuszczenia były słuszne - zastosowanie entropii krzyżowej przynosi dużo lepsze rezultaty niż zastosowanie funkcji MSE. W przypadku dokładności predykcji na zbiorze treningowym różnica między nimi wynosi ponad 30% accuracy.

Jak wspomnieliśmy, takie rezultaty mogą wynikać z kilku przyczyn. Po pierwsze, cross entropy jest lepszą miarą niż MSE do klasyfikacji, ponieważ granica decyzyjna w zadaniu klasyfikacyjnym jest duża (w porównaniu z regresją). MSE nie "penalizuje" wystarczająco błędnych klasyfikacji, ale za to jest właściwą funkcją błędu w przypadku regresji, gdzie odległość między dwiema wartościami, które można przewidzieć, jest niewielka.

Po drugie, z probabilistycznego punktu widzenia, entropia krzyżowa to oczywista funkcja kosztu, w przypadku, jeśli w warstwie wyjściowej sieci występuje nieliniowość sigmoidalna (lub softmax) i chcemy zmaksymalizować prawdopodobieństwo prawidłowego sklasyfikowania danych wejściowych. Do klasyfikacji entropia krzyżowa wydaje się być bardziej odpowiednia niż MSE - podstawowe założenia po prostu mają więcej sensu dla tej "konfiguracji". Mówiąc to, mamy na myśli, że oczywiście można wyszkolić klasyfikator z funkcją błędu MSE i prawdopodobnie będzie działał dobrze, jednak jak widzimy - nie aż tak dobrze, jak ma to miejsce w przypadku cross entropy. W celu poprawy wyników MSE lepszym pomysłem na funkcję aktywacji byłaby prawdopodobnie funkcja liniowa.

Konkludując - cross entropy świetnie sprawdza się w rozwiązywaniu problemów klasyfikacyjnych -czyli właśnie tak jak ma to miejsce w przypadku naszego zagadnienia. W przypadku problemów regresyjnych lepszym wyborem powinno być MSE.

Analiza dokładności sieci w zależności od liczby iteracji

Epoka - nazywana w naszej sieci iteracją oznacza uczenie sieci neuronowej ze wszystkimi danymi treningowymi przez jeden cykl. W epoce używamy wszystkie dane dokładnie raz. Forward propagation i backpropagation liczone są łącznie jako jedno "przejście":

Epoka składa się z co najmniej jednego batcha, gdzie używamy części zbioru danych do uczenia sieci neuronowej. Przejście przez przykłady treningowe w batchu nazywamy iteracją.

Czasami nazwa epoka jest używana wymiennie z "iteracją". Aby jednak w pełni wyjaśnić pojęcia, rozważymy prosty przykład, w którym mamy 1000 punktów danych, jak pokazano na poniższym rysunku:

Jeśli rozmiar batcha wynosi 1000, możemy zakończyć epokę pojedynczą iteracją. Podobnie, jeśli rozmiar batcha wynosi 500, epoka wymaga dwóch iteracji. Tak więc, jeśli rozmiar batcha wynosi 100, ukończenie epoki zajmuje 10 iteracji.

Wracając jednak do meritum - podjęcie decyzji o architekturze sieci neuronowej to duży krok w budowie modelu. Musimy w końcu wytrenować model i dostroić po drodze więcej hiperparametrów.

Na etapie treningu dążymy do zminimalizowania błędu, a także upewnienia się, że model dobrze uogólnia nowe dane. Kompromis między obciążeniem a wariancją jest nadal potencjalną pułapką, której chcemy uniknąć, podobnie jak w przypadku innych nadzorowanych algorytmów uczenia maszynowego.

Z tzw. overfittingiem, czyli nadmiernym dopasowaniem (duża wariancja) mamy do czynienia, gdy model idealnie pasuje do przykładów treningowych, ale ma ograniczone zdolności do uogólnienia. Z drugiej strony, jeśli mówi się, że model jest niedopasowany (duże obciążenie estymatora), jeśli nie nauczył się danych wystarczająco dobrze:

Oczekuje się, że dobry model uchwyci podstawową strukturę danych. Innymi słowy, nie dochodzi do overfittingu jak również underfittingu.

Budując model sieci neuronowej, ustawiamy parametr liczby epok przed rozpoczęciem uczenia. Jednak początkowo nie możemy wiedzieć, ile epok jest korzystnych dla modelu. W zależności od architektury sieci neuronowej i zestawu danych musimy zdecydować, kiedy wagi sieci neuronowych są zbieżne.

W przypadku modeli sieci neuronowych często bada się wykresy krzywych uczenia, aby zdecydować o zbieżności modelu. Ogólnie rzecz biorąc, tworzymy wykresy błędu lub dokładności w danej epoce. Podczas treningu spodziewamy się, że strata zmniejszy się, a dokładność wzrośnie wraz ze wzrostem liczby epok. Oczekujemy jednak, że po pewnym czasie zarówno błąd, jak i dokładność ustabilizują się.

Zwykle zaleca się podzielenie zbioru danych na zbiory uczące i walidacyjne - tak jak zrobiliśmy też w naszym przypadku. W ten sposób możemy stworzyć wykresy krzywej uczenia dla różnych zestawów. Te wykresy pomagają nam zdiagnozować, czy model jest nadmiernie wyuczony, niedouczony lub pasuje do zestawu uczącego.

W zależności od architektury i dostępnych danych, liczbę epok możemy traktować jako hiperparametr.

Poniżej przedstawione zostały właśnie takie wykresy:

Jak już wspomnieliśmy, wagi sieci neuronowych są aktualizowane iteracyjnie, ponieważ jest to algorytm oparty na gradient descent. Pojedyncza epoka treningu nie wystarczy i prowadzi do niedopasowania. Biorąc pod uwagę złożoność rzeczywistych problemów, wytrenowanie sieci neuronowej może zająć wiele epok. W rezultacie spodziewamy się, że wykresy krzywej uczenia będą coraz lepsze, aż do osiągnięcia konwergencji. Następnie, jeśli będziemy dalej trenować model, będzie on nadmiernie dopasowany, a błędy walidacji zaczną się zwiększać.

Patrząc na powyższe wykresy, możemy zauważyć, że dla naszej sieci najgorsze rezultaty odnosi liczba epok równa 8. Jak pisaliśmy, wynika to z faktu, że jeśli ustawimy zbyt małą liczbę epok, to uczenie zatrzyma się jeszcze przed zbieżnością modelu.

Najlepsze rezultaty zaś są widoczne dla liczby epok równej 20 oraz 25. Dla naszego modelu jest to wystarczająca ilość, by spełnić kompromis między nauką a generalizacją.

Przegląd literatury

Co ciekawe, warto zauważyć, że nasze spostrzeżenia dotyczące accuracy modelu w zależności od jego parametrów co do meritum pokrywają się z obserwacjami zawartymi w artykule "On Machine Learning Models for Heart Disease Diagnosis" [Chu-Hsing Lin, Po-Kai Yang, Yu-Chiao Lin, Pin-Kuei Fu], choć nieco się od niego różnią. We wspomnianej pracy pokazano porównanie modelu NN wśród różnych wariantów liczby warstw ukrytych: od 1 do 3, a także w różnym zakresie liczby neuronów na warstwę: mniej niż 20; od 21 do 100, ponad 100. Wyniki zostały przedstawione w poniższej tabeli.

Accuracy of NN on variant numbers of hidden layers and neurons in each layer

NN on variant numbers of hidden layers and neurons in each layer

Mówiąc, że nasze spostrzeżenia zgadzają się "co do meritum", mamy na myśli, że więcej nie zawsze znaczy lepiej. W przypadku naszych specyficznych parametrów korzystniejsze wyniki odnieśliśmy jednak z użyciem ilości neuronów mniejszej od 100 jednak większej niż 20 - konkretnie 64 oraz liczby warstw równej 3.

Inne powiązane artykuły naukowe

Istnieje wiele wcześniejszych prac nad klasyfikacją chorób serca. W jednym z wczesnych badań (1989 r.) accuracy 77,0% zostało uzyskane przy użyciu algorytmu regresji logistycznej ("International application of a new probability algorithm for the diagnosis of coronary artery disease" [R. Detrano et al]). W tym samym roku, w innym badaniu, w którym autor wykorzystał koncepcyjny system grupowania CLASSIT, wynik trafności wyniósł 78,9% ("Models of incremental concept formation" [J. H. Gennari, P. Langley and D. Fisher, ]. Zaobserwowano również, że accuracy na poziomie 81,11% i 81,48% można było osiągnąć, stosując odpowiednio algorytmy C4.5 i Naïve Bayes ("Machine learning techniques for medical analysis. School of Information Technology and Electrical Engineering" [N. Cheung] 2001r.).

Jak zatem prezentują się obecnie najlepsze wyniki dla naszego zagadnienia? Poniżej przedstawimy odpowiedź to pytanie.

Jak moglibyśmy poprawić zdolności predykcyjne naszego modelu - analiza publikacji naukowej

W artykule "Classification models for heart disease prediction using feature selection and PCA" [Anna Karen Gárate-Escamila, AmirHajjam El Hassania, Emmanuel Andrèsbc] zaproponowano model oparty na metodzie chi-kwadrat i PCA do problemu wykrywania chorób serca. Wyniki eksperymentalne wykazały, że PCA zapewnia lepszą prognozę dotyczącą problemu klasyfikacji wielowymiarowej poprzez dobór cech metodą chi-kwadrat. Chi-kwadrat klasyfikuje zmienne niezależne najbardziej zgodne z etykietą. We wspomnianej pracy uzyskane wyniki (za pomocą metody CHI-PCA-RF) wynosiły 99.0% accuracy, 100.0% precyzji.

Biorąc pod uwagę te rezultaty, uważamy, że roztropnym byłoby rozważenie metody analizy głównych składowych przy budowie kolejnych sieci neuronowych.

Analizę porównawczą otrzymanych przez nas wyników z innymi metodami

Jak wspomnieliśmy powyżej, zagadnienie wykrywania chorób serca może być - jak większość problemów - rozwiązane na wiele sposobów.

Jednym z nich może być użycie CNN jak w artykule "On Machine Learning Models for Heart Disease Diagnosis" [Chu-Hsing Lin, Po-Kai Yang, Yu-Chiao Lin, Pin-Kuei Fu].

Jak możemy dowiedzieć się z publikacji, konwolucyjne sieci neuronowe (CNN) mają inną architekturę w porównaniu ze zwykłymi sieciami NN w kilku aspektach. Neurony w jednej warstwie nie łączą się ze wszystkimi neuronami w warstwie poprzedniej. Ponadto CNN mają dwie części składowe: ekstraktor cech i klasyfikator. W części ekstraktora cech sieć wykonuje serię operacji splotowych (przez warstwy splotowe) i ograniczeń (przez warstwy łączące), za pomocą których cechy są wykrywane.

Oto jak wyglądają otrzymane przez nich rezultaty:

Jak zauważają sami autorzy dokładność modelu CNN nie jest jednak tak wysoka, jak w konwencjonalnym modelu NN, ale w prawie wszystkich przypadkach algorytm zachowywał się stabilnie pod względem wydajności. W zagadnieniu wykrywania chorób serca CNN wykazał gorsze wyniki niż sieć z warstwami w pełni połączonymi, a jednym z możliwych powodów, jest rozmiar zbioru danych, wynoszących 303 rekordy. Jak wspominają autorzy publikacji, w przyszłości należałoby sprawdzić poprawność modelu dla większej ilości danych.

Warto w tym miejscu również dodać, że CNN jest algorytmem używanym najczęściej do klasyfikacji obrazów. Wyniki accuracy na poziomie 82% są zatem dla nas zaskakująco dobre, gdyż nie spodziewaliśmy się, że CNN może się sprawdzić w naszym zagadnieniu. Można zatem oczekiwać, że autorzy publikacji posługiwali się CNN, w których filtry były jednowymiarowe.

Poniżej przedstawiamy również wyniki uzyskane za pomocą innych metod w artykule "Neural network diagnosis of heart disease" [Ebenezer O. Olaniyi; Oyebade K. Oyedotun; Abdulkader Helwan; Khashman Adnan]

Niestety trudno jest nam odnieść się do metod działania wszystkich algorytmów oraz wyników ich dokładności, gdyż nie mamy jeszcze wystarczającego doświadczenia w tej dziedzinie. Niemniej jednak stwierdziliśmy, że warto wspomnieć o skuteczności innych metod, jak i również pokusiliśmy się o samodzielną implementację paru z nich.

Poniżej prezentujemy zatem metody, z którymi spotkaliśmy się do tej pory i uznaliśmy, że również mogą posłużyć do rozwiązania naszego problemu. Przedstawione poniżej modele zostały stworzone przez nas samodzielnie w programie Statistica:

Krokowa wsteczna:

Najlepszy podzbiór

Krokowa postępująca

Probit – wszystkie efekty

Probit – krokowa postępująca

Probit – krokowa wsteczna

Probit – najlepszy podzbiór

W modelach logistycznym oraz probitowym w naszym zbiorze zmiennych zmienne: age, trestbps, chol, thalach oraz oldpeak były traktowane jako zmienne ilościowe, natomiast pozostałe zmienne jako jakościowe (w koniecznych przypadkach kodowane za pomocą kilku dummies).

Zwróćmy uwagę, że są to modele oszacowane dla całego zbioru danych. Gdybyśmy przeprowadzili tu podział na zbiór uczący i testowy, rozsądnie byłoby przyjąć, że wyniki mogłyby się pogorszyć lub nie pogorszyć, ale na pewno nie poprawić (pomijając przypadek gdybyśmy wiele godzin pracy włożyli w odpowiednią preselekcię zbioru uczącego i testowego, co jednak nie poprawiłoby ogólnych zdolności predykcyjnych modelu i byłoby czystym artefaktem obliczeniowym).

Oszacowaliśmy kilka modeli logistycznych (w oparciu o jedno ze standardowych kryteriów: statystyki punktowej, funkcji największej wiarygodności oraz kryterium Akaike’a) oraz kilka probitowych w klasycznych scenariuszach uwzględniania wszystkich zmiennych, wyboru najlepszego podzbioru zmiennych oraz regresji krokowej postępującej i wstecznej.

Wyniki tych estymacji były podobne. Regularnie nieco trafniej były przewidywane przypadki śmierci (=1) niż przypadki przeżycia (=0), ale tego można było się spodziewać ze względu na niedoskonałe zrównoważenie próby w zmiennej objaśnianej (54% jedynek).

Można by to próbować korygować doprowadzając do zrównania się trafności dla obu podklas, nie robiliśmy jednak tego, bo wyniki były i tak słabsze niż rezultaty osiągane przez naszą sieć.

Ogólne spostrzeżenia są takie, że w tej klasie modeli modele logitowe były nieco bardziej trafne wobec probitowych. Trafność kształtowała się w przedziale od 83% do 93%. Naturalnie najbardziej trafne zawsze były modele z pełnym zestawem zmiennych objaśniających (także tych nieistotnych statystycznie na umownym poziomie istotności 0,05).

Realnej trafności na całkowicie zewnętrznych danych intuicyjnie można by oczekiwać zapewne gdzieś w okolicach 80%-85%, czyli nieco poniżej dopasowani modelu tylko ze zmiennymi istotnymi statystycznie.